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Hysteresis in the Random Field Ising Model and Bootstrap Percolation 
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We study hysteresis in the random-field Ising model with an asymmetric distribution of quenched fields, in 
the limit of low disorder in two and three dimensions. We relate the spin flip process to bootstrap percola- 
tion, and show that the characteristic length for self-averaging L* increases as exp(exp(J/A)) in 2d, and as 
exp(exp(exp(J/A))) in 3d, for disorder strength A much less than the exchange coupling J. For system size 
1< L < I*, the coercive field h collr varies as 2 J — A In In L for the square lattice, and as 2 J — A In In In L 
on the cubic lattice. Its limiting value is for L — > oo both for square and cubic lattices. For lattices with 
coordination number 3, the limiting magnetization shows no jump, and h coer tends to J. 
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In recent years, there has been a lot of interest in the 
study of hysteresis in magnetic systems, both theoretically 
and in experiments [^J. Hysteresis in the random-field Ising 
model (RFIM) model was first discussed by Sethna et al [|3|], 
who proposed it as a model of return point memory, and of 
Barkhausen noise Q. Sethna et al solved the model in the 
mean-field limit, and showed that if the strength A of the 
quenched random field is large, the average magnetization 
per site is a continuous function of the external field, but for 
small A, it shows a discontinuous jump as the external field 
is increased. Interestingly, the non-equilibrium hysteresis re- 
sponse in the RFIM can be determined exactly on a Bethe lat- 
tice though the corresponding equilibrium problem has 
not been solved so far, even in zero field. These calculations 
have been extended to determining the distribution of sizes of 
the Barkhausen jumps [Q], and the calculation of minor hys- 
teresis loops 

In this paper, we study the low disorder limit of the hys- 
teresis loop in the RFIM on periodic lattices in two and three 
dimensions. We find that the behavior of hysteresis loops de- 
pends nontrivially on the coordination number z. For z = 3, 
for continuous unbounded distributions of random fields, the 
hysteresis loops show no jump discontinuity of magnetization 
even in the limit of small disorder, but for higher z they do. 
This is exactly as found in the exact solution on the Bethe 
lattice [§]. 

The analytical treatment of self-consistent equations on the 
Bethe lattice is immediately generalized to asymmetrical case. 
However, we find that behavior of hysteresis loops in eu- 
clidean lattices can be quite different from that on the Bethe 
lattice, for asymmetrical distributions. On hypercubical lat- 
tices in d dimensions, there is an instability related to boot- 
strap percolation, that is absent on the Bethe lattice. This re- 
duces the value of the coercive field h cocr away from the Bethe 
lattice value 0( J) to zero, where J is the exchange coupling. 
We note that the limit A — > is somewhat subtle, as the sys- 
tem size L* required for self-averaging diverges very fast for 
small A, and the finite-size corrections to the thermodynamic 
limit tend to zero very slowly. 

In the RFIM, the Ising spins {si} with nearest neighbor fer- 
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romagnetic interaction J are coupled to the on-site quenched 
random magnetic field hi and the external field h. The Hamil- 
tonian of the system is given by 



H = -J 2J s i s j - hjSj - h , 



(1) 



<»,i> 



We assume that {hi} are quenched independent identically 
distributed random variables with the probability that the 
value of the random field at site i lies between hi and hi + dhi 
being cf)(hi) dhi. 

The system evolves under the zero-temperature Glauber 
single-spin-flip dynamics [|Io|: a spin-flip is allowed only if 
the process lowers energy. We assume that the rate of spin- 
flips is much larger than the rate at which h is changed, so 
that all flippable spins may be said to relax instantly, and any 
spin Si always remains parallel to the net local field li at the 
site: 



sign( 



i) = sign(J^ Sj + hi + h). 



(2) 



Under this dynamics, for ferromagnetic coupling (J > 0), 
if we start with any stable configuration, and then increase the 
external field and allow the system to relax, the final stable 
configuration reached is independent of the order in which the 
unstable spins are flipped. Also, in the relaxation process no 
spin flips more than once. 

For a given distribution <fi(hi), we define p m (h) with < 
to < z as the conditional probability that the local field at 
any site i will be large enough so that it will flip up, if to of its 
neighbors are up, when the uniform external field is h. Clearly 



Pm{h) = 



(z — 2m)J—h 



<t>(hi) dhi. 



(3) 



Clearly, for any given value of h, the magnetization depends 
on the distribution <j>(hi) only through p m (h). 

Historically, RFIM was first studied in the context of pos- 
sible destruction of long range order by arbitrarily weak 
quenched disorder in equilibrium systems. Accordingly the 
distribution of random field was assumed to be symmetrical. 
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However, in hysteresis problem, the symmetry between up 
and down spins state is already broken by the specially pre- 
pared initial state ( all down in our case), and the symmetry 
of the distribution plays no special role. In the following, we 
shall assume that the distribution has a asymmetrical shape, 
given by 

4>{hi) = ^exp(-hi/A)6(hi); (4) 

where 8 is the step function. The mean value of hi can be 
made zero by a shift in the value of the external uniform field. 
Our treatment is easily extended to other continuous unimodal 
distributions. The exact form of 4>{x) is not important, and 
other forms like cxp(— x — e~ x ) which fall sharply for nega- 
tive x have the same behavior. 

Consider first the case of the two-dimensional hexagonal 
lattice with z = 3. For periodic boundary conditions (PBC), 
if A = 0, starting with a configuration with all spins down, 
clearly one has h cocr = 3 J. For A ^ 0, the site with the 
largest local field flips first, and then if h > J, pi (h) = 1, this 
causes neighbors of the flipped spin to flip, and their neigh- 
bors, and so on. Thus, so long as there is at least one flipped 
spin, all other spins also flip, and the magnetization is 1. The 
largest local field in a system of L 2 spins is of order 2 A In L. 
Once this spin turns up, other spin will flip also up, causing 
a jump in magnetization from a value w —1 to a value +1 in 
each sample. Hence the coercive field, (the value of h where 
magnetization changes sign) to lowest order in A, is given by 

h cocr = 3 J - 2A In L, for 1 < In L < J /A. (5) 

Sample to sample fluctuations in the position of the jump are 
of order A. On averaging over disorder, the magnetization 
will become a smooth function of h, with the width of the 
transition region being of order A. 

For a fixed A <C J, if L is increased to a value near 
exp(J/A) = L* hex , h coer decreases to a value near J. For 
h as J, pi(h) is no longer nearly 1, but po(h) ~ 0,p2(h) ~ 
Ps(h) ~ 1. The value of magnetization depends only on 
Pi{h), which is a function of h = (h— J)/ A. As h increased 
from — oo, Pi(h) increases continuously from to 1. 

In Fig. [j], curve A shows the result of a simulation on the 
hexagonal lattice with L = 4096, and PBC. To avoid the prob- 
lem of probability of nucleation being very small for h near J, 
we made the local field at a small fraction of randomly cho- 
sen sites very large, so that these spins are up at any h. The 
number of such spins is of order L, so that their effect on 
the average magnetization is negligible. Introduction of these 
"nucleation centers" makes L* sw 0(y/~L) ( the average sep- 
aration between centers), and h cocr drops to a value near J, 
so that, we can study the large L limit with available comput- 
ers. For L > L^ ex , the behavior of hysteresis loops becomes 
independent of L. 

We see that magnetization no longer undergoes a single 
large jump, but many small jumps. In the figure, we also show 
the plot of magnetization when the random field at each site 



is decreased by a factor 10. This changes the value A from 
0.1J to 0.01J. However, plotted as a function of h, the mag- 
netization for these two different values (for small A) fall on 
top of each other for the same realization of disorder (except 
for the overall scale A). Thus we can decrease A further to 
arbitrarily small values, and the limit of A — > is straight- 
forward for each realization of disorder. Then, averaging over 
disorder, for a fixed A, we see that h coer tends to the value 
J as A tends to 0. Also, we see that there is no macroscopic 
jump-discontinuity for any non-zero A. 

We also show in Fig. [j] [curve B], the results of simulation 
of a 3-dimensional lattice with z = 3 of size 256 3 with PBC. 
The behavior is qualitatively same as that in two dimensions. 
The value of h coer — J in the limit A — > is same for sym- 
metrical distribution, and also is the same as predicted by the 
Bethe approximation. 

On the square lattice also, the value of h cocr is determined 
by the need to create a nucleation event. Arguing as before, 
we see that h coer to lowest order in A is given by h coer ~ 
4 J - 2A In L, for l<lnL< J / A. Adding a small number 
of nucleation sites suppresses this slow transient, and lowers 
/i cier from 4 J to a value near 2 J. However, in this case, even 
after adding the nucleation centers, the system shows a large 
single jump in magnetization, indicating the existence of an- 
other instability. We observed in the simulation that at low A, 
as h is increased, the domains of up spins grow in rectangu- 
lar clusters [see Fig. ^ and at a critical value of h coer , one of 
them suddenly fills the entire lattice. This value h coer fluctu- 
ates a bit from sample to sample. In Fig. |3] we have plotted 
the distribution of the scaled variable h c — {h coer — 2J)/A 
for different system sizes L, for A = 0.001J. The number of 
different realizations varies from 10 4 (for the largest L) to 10 5 
(for the smallest L). Note that the distribution shifts to the left 
with the increasing system size, and becomes narrower. 

This instability can be understood as follows: on a square 
lattice, for the asymmetric distribution [Eq. for h > 0, 
Pm = 1 for to > 2, and any spins with more than one up- 
neighbors flips up. Therefore, stable clusters of up spins are 
rectangular in shape. The growth of domains of up spins is 
same as in the bootstrap percolation process BP m with m = 2 
|[TT|— [l~3|] . In the process BP m , the initial configuration is pre- 
pared by occupying lattice sites independently with a proba- 
bility p and the resulting configuration is evolved by the rules: 
the occupied sites remain occupied forever, while an unoc- 
cupied site having at least to occupied neighbors, becomes 
occupied. For m — 2, on a square lattice, in the final config- 
uration, the sites which are occupied form disjoint rectangles, 
like the cluster of up-spins in Fig. ||. It has been proved that 
in the thermodynamic limit of large L, for any initial concen- 
tration p > 0, in the final configuration all sites are occupied 
with probability 1 JT^]. 

Now consider a rectangular cluster of up spins, of length 
I and width to. Let P(l, to) be the probability that, if this 
rectangle is put in a randomly prepared background of den- 
sity pi(h), this rectangle will grow by the BP2 process to fill 
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the entire space. The probability that the random fields at any 
sites neighboring this rectangle will be large enough to cause 
it to flip up is pi(h). The probability that there is at least one 
such site along each of two adjacent sides of length I and m of 
the rectangle is (1 — <r)(l — q m ), where q = l—pi(h). Once 
these spins flip up, this induces all the other spins along the 
boundary side to flip up and the size of the rectangle grows to 
(/ + 1) x (m + 1). Therefore 



P{l,m) > (l-q l )(l-q m )P(l + l,m + l). 



(6) 



Thus the probability of occurrence of a nucleation which fi- 
nally grows to fill the entire lattice is 



Pnuc>Po(h)Y[(l 

3=1 



^ 2 



(7) 



The right hand side can be shown to vary as 
Po(^) ex P (— zp x (h) ) for small pi (ft.). The condition to deter- 
mine h coer is that for this value of ft, P nU c becomes of order 
1 /J 2 , so that we get 



Po{h coer ) exp 



3pi(fc„„.r) 



1 

Z 2 ' 



(8) 



This equation can be solved for h coer for any given L. For the 
distribution given by Eq. (0), this becomes 



exp 



ft -4J 



exp 



2vr 2 /-ft 
— — exp 



2 J 



1 

'' Z 2 ' 
(9) 



It is easy to see from this equation that for 
1 -C InL <SC J/A, the leading i-dependence of h cocr , 
to lowest order in A is given by 



h 



4 J — 2AlnL; 



and for J/A < InL < exp(2J/A), 

ft cocr «2J- Aln[-^(lnL- J/A)]. 



(10) 



(11) 



This agrees with our observation that the scaled critical field 
h c shifts to the left with increasing system size. 

To test the validity of Eq. (p) in simulations, we put po(h) = 
0.005 independent of ft. Eq. (8h then simplifies to 



pi(K o J 



6lnL 



(12) 



In Fig. |], we have plotted p\ for the mean h coer from Fig. || 
versus 1 / In L. The graph is approximately a straight line, 
which agrees with Eq. (^2j). The slope of the line is less than 
in Eq. (p"2|), which only gives an upper bound to /i eocr . 

If h > 0, we will have p2 — 1, and bootstrapping ensures 
that so long as po > 0, we will have all spins up in the limit 
of large L. This implies that h coer = in this limit. 



If there are sites with large negative quenched fields, the 
bootstrap growth stops at such sites. Hence the bootstrap in- 
stability cannot be seen for symmetric distributions. Even if 
the quenched fields are only positive, the instability does not 
occur on lattices with z = 3. On such lattices, if the un- 
occupied sites percolate, there are infinitely extended lines of 
unoccupied sites in the lattice. These cannot not become occu- 
pied by bootstrapping under BP2. Thus the critical threshold 
for BP2 on such lattices is not 0. 

The above analysis is easily extended to higher dimensions. 
In d = 3, if ft > 0, then p m (h) = 1 for m > 3, therefore the 
spin flip process is similar to the spanning process of three 
dimensional BP3 jl4j]. But in this case, it is known that for 
any initial non-zero density, in the thermodynamical limit, the 
final configuration has all sites occupied with probability 1. 
The clusters of up-spins grow as cuboids, and at each surface 
of the cluster, the nucleation process is similar to that in two 
dimension. Let e be the probability that, a nucleation occurs 
at a given point of a surface of the clusters of up spins which 
sweeps the entire two dimensional plane at ft. 



Pi (ft) exp 



3p2(h)J 



(13) 



The probability that, there exist at least one nucleation which 
sweeps the entire plane of size / x /, is 1 — (1 — e) . Therefore, 
the probability P nuc , that a nucleation sweeps the entire three 
dimensional lattice at h satisfies 



oc 

Pnuo>Po(h)l[[l-(l-e) 1 



1=1 



(14) 



For small e, the infinite product can be shown to vary as 
exp(-A/ v / e), with A = §^(3/2). 



h coer is determined by the condition that P n 
the order 1/J 3 : 



Po(ft coer )exp 



A 



: exp 



must be of 

a 1/J 3 . 
(15) 



The leading L-dependence of ft cocr is different in dif- 
ferent ranges of h coer , depending on whether the strongest 
dependence of the left-hand side comes from variation of 
p {h),pi(h) or p 2 {h). We find that, h coar w 6 J - 3 A InL, 
for 4 J < ft cocr < 6 J. It is « 4 J - 2Aln(lni - 2J/3A), 
for 2 J < h°"„ < 4 J; and « 2 J - Alnln(ln£ - 2J/3A), 
for A <C h rocT < 2 J. It is straightforward to determine the 
corresponding ranges of L for the validity of these equations. 

In the limit L ^> L* ub = exp(exp(exp(2J/A))), the loop 
becomes independent of L, with h cger — > 0. We have also ver- 
ified the existence of jump in numerical simulation for z = 4 
(diamond lattice) in three dimensions. 

In brief, we have shown that the hysteresis loops on lat- 
tices with coordination number three are qualitatively differ- 
ent from those with z > 3. For the square and cubic lattices, 
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h eoer decreases to A very slowly for large L. In general, it is 
true for lattices where the corresponding bootstrap percolation 
problem has an instability. 

We thank M. Barma and N. Trivedi for critically reading 
the manuscript. 
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FIG. 2. A snapshot of the up-spins just before the jump 
(h = 1.998243J). The lattice size is 200 x 200 and A = 0.001 J. 
Initial configuration is prepared with 0.05% up-spins. 
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FIG. 3. Distribution of the scaled coercive field on a square lattice 
for different lattice size L 2 . 




